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CONVERSION  FACTORS,  NON-SI  TO  SI  (METRIC) 
UNITS  OF  MEASUREMENT 


Non-SI  units  of  measurement  used  in  this  report  can  be  converted  to  SI  (metric) 
units  as  follows: 


Multiol 


To  Obtain 


inches 

2.54 

centimetres 

kilofeet 

0.3048 

kilometres 

kilotons 

(nuclear  equivalent  of  TNT) 

4.184 

tera joules 

pounds  (mass) 

0.45359237 

kilograms 

pounds  per  square  inch 

6.894757 

kilopascals 

pounds  per  square  inch 
per  second 

6.894757 

kilopascals/second 

tons 

(nuclear  equivalent  of  TNT) 

4.184 

giga joules 

RESPONSE  OF  AN  ELASTIC-PLASTIC  STRUCTURE  SUBJECTED  TO 
A  SIMULATED  NUCLEAR  BURST 


CHAPTER  1 
INTRODUCTION 

Very  few  structures  are  designed  to  remain  elastic  under  highly  impulsive 
loads.  Systems  that  are  designed  to  behave  elastically  under  such  loads  are 
generally  uneconomical  and  tend  to  have  structural  characteristics  that  at¬ 
tract  greater  dynamic  forces.  Consequently,  some  plastic  behavior  is  desir¬ 
able  for  the  largest  dynamic  loads  that  are  anticipated  to  occur. 

The  stress  distribution  in  a  structure  that  undergoes  a  plastic  response 
is  complex  and  difficult  to  predict.  Hence,  the  design  and  failure  criteria 
for  such  structures  are  generally  governed  by  limitations  on  respohse,  e.g., 
deflections  and  rotations,  rather  than  by  limitations  on  stresses. 

This  paper  presents  the  mathematical  model  and  the  numerical  solution  for 
the  response  of  structures  that  behave  as  elastic-perfectly  plastic  systems 
when  subjected  to  monotonically  decaying  loads.  The  loading  definition  is  in¬ 
tended  to  simulate  the  pressure-time  history  of  a  nuclear  pulse  with  zero  rise 
time  and  no  negative  phase. 


CHAPTER  2 
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ANALYTICAL  PROCEDURE 

2.1  DEFINITIONS  OF  DUCTILITY 

AND  RESISTANCE 

A  Measure  of  the  total  elastic  and  plastic  deformation  in  a  structure  is 

its  ductility.  In  this  paper,  ductility  will  be  defined  as  the  ratio  of  the 

maximum  displacement  under  a  prescribed  load  to  peak  elastic  displacement, 

i.e.,  x  /x  1  in  Figure  2.1.  The  maximum  ductility  that  a  structure  can  sus- 
m  y 

tain  without  collapse  is  controlled  by  its  material  and  geometric  characteris¬ 
tics  and  is  an  indicator  of  the  amount  of  energy  that  can  be  absorbed  prior  to 
failure. 

The  relationship  between  static  load  and  deflection,  i.e.,  the  resistance 
function,  for  many  structural  elements  may  be  idealized  as  the  elastic- 
perfectly  plastic  curve  shown  in  Figure  2.1.  For  monotonic  dynamic  loads,  it 
is  generally  assumed  that  the  resistance  function  maintains  the  same  charac¬ 
teristic  form  as  for  the  static  situation.  However,  the  magnitudes  may  be 
increased  to  represent  the  enhanced  strength  of  the  material  under  the  high 
strain  rates  associated  with  dynamic  loads. 

2.2  DESCRIPTION  OF  LOAD 

In  the  past,  nuclear  pressure  histories  have  often  been  represented  as 
simple  triangular  functions  or  as  superimposed  combinations  of  simple  tri¬ 
angles  (References  1  and  2).  Analyses  made  with  such  loading  definitions  have 
often  proved  to  be  reliable,  particularly  when  the  loads  are  most  impulsive  or 
when  the  duration  of  the  load  is  large  compared  to  the  natural  period  of  the 

structure.  However,  when  attempting  to  represent  some  nuclear  weapon  load¬ 

ings,  the  use  of  an  "equivalent"  triangular  function  introduces  significant 
error.  The  selection  of  an  appropriate  triangle  to  represent  a  particular  nu¬ 
clear  pulse  adds  one  more  degree  of  uncertainty  to  the  already  subjective  dy¬ 

namic  analysis. 

Analytic  approximations  to  actual  nuclear  burst  overpressures  have 
recently  been  developed.  The  most  recent  expression  for  the  overpressures  can 
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For  convenience,  symbols  and  unusual  abbreviations  are  listed  and  defined  in 
the  Notation  (Appendix  C). 


be  found  in  Reference  3  with  subsequent  modifications  by  memoranda  from  the 
authors.  The  Speicher-Brode  fit  relates  height  of  burst,  range,  time  of  ar- 
rival,  positive  phase  duration,  and  overpressure  for  a  1-kiloton  (KT)  weapon. 
Calculation  of  the  quantities  listed  above,  as  well  as  positive  phase  impulse, 
for  a  1  megaton  (NT)  surface  burst  weapon  are  included  in  Table  2.1  and  Figure 
2.2.  Cube-root  scaling  can  be  used  to  relate  the  quantities  of  interest  for 
different  weapon  yields. 

2.3  STRUCTURAL  IDEALIZATION 

Many  structures  may  be  idealized  as  single-degree-of-freedom  (SDOF) 
spring-mass  models  in  order  to  determine  their  fundamental  response.  More 
complex  structures  may  be  idealized  as  a  combination  of  spring~mass  models 
that  are  superimposed  by  modal  superpositioning.  The  idealized  SDOF  model 
used  in  this  paper,  as  well  as  the  resistance  and  loading  idealizations,  are 
shown  in  Figure  2.1. 

The  differential  equation  of  motion  for  the  undamped  SDOF  model  is  given 
below. 

mx  +  kx  =  F(t)  (2.1) 

where 

2  2 

x  =  d  y/dt  =  acceleration 

and  kx  =  r  ,  if  plastic  response  has  occurred 
in 

Viscous  damping  is  generally  not  included  in  the  analysis  of  structures 
that  are  subjected  to  short-duration  monotonic  loads,  particularly  when  only 
the  peak  response  is  of  primary  interest.  Damping  at  less  than  5  or  6  percent 
of  critical  will  have  only  a  minor  effect  on  peak  response  for  nonoscillating 
loads . 

Equation  2.1  can  be  nondimensionalized  by  a  transformation  of  the  inde¬ 
pendent  variables  to  yield  the  following  equations  from  Reference  2. 

i  Fi 

— ^  rj  +  n  =  —  f (4)  ,  if  elastic,  (2.2) 

Un  rm 

^A  table  of  factors  for  converting  U.  S.  customary  units  of  measurement  to 
metric  (SI)  units  is  presented  on  page  4. 


and 


1  *  1 

~2  rj  +  1  =  ^  f(4)  »  if  plastic, 


(2.3) 


where 


n(4)  =  x(t)/xy 

rj(0  =  T2x(t)/xy 

1  =  t/T 

xy  =  displacement  at  first  yield 


T  =  2n 


4T 


=  natural  period 


It  should  be  noted  that  the  specified  forcing  function  in  Equations  2.2 
and  2.3  implies  that  pressures  may  be  normalized  with  respect  to  a  single  time 
parameter.  However,  the  Speicher-Brode  loading  definition  is  a  function  of 
time  and  of  peak  overpressure  and  may  not  be  normalized  by  a  single  parameter 
as  illustrated  in  Figure  2.3.  Thus,  a  separate  response  chart  must  be  devel¬ 
oped  for  each  value  of  peak  overpressure,  P  . 

s 

2.4  NUMERICAL  SOLUTION 

There  are  several  numerical  methods  for  solving  the  nondimensionalized 
equations  of  motion.  A  method  that  is  frequently  applied  to  problems  in  struc¬ 
tural  dynamics  is  the  Newmark  p  method  (Reference  4).  This  method  was  se¬ 
lected  for  numerical  analysis  because  of  its  flexibility  and  established  cri¬ 
teria  for  stability  and  convergence. 

From  previously  defined  relationships,  we  have 

x(t)  =  xyn(£) 
x(t)  =  xyiK4)/T 

x(t)  =  xyn(U/T2 


t  =  - 
s  T 


For  our  initial  conditions 


n0  =  0 
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1*1  =  0 
'o 

The  initial  value  of  the  acceleration  term  is 

fi  =  4n2F,/r  _ _ 

'o  r  m  9 

The  Newmark  P  equations  may  be  written  in  nondimensionalized  form  as 

indicated  below: 


Vi  *  \  *  (r)  \  -  « (r)\  ♦ » (r)2  «-5> 


V,  *  <ii  =  w  (r)^  -  «•« 

By  making  an  assumption  for  the  acceleration  at  the  (i+1)  time  step,  ap¬ 
proximations  of  the  corresponding  displacement  and  velocity  terms  can  be  made 
with  the  equations  above.  The  resulting  displacement  can  then  be  inserted  in¬ 
to  the  original  differential  equation  and  tested  for  convergence.  The  process 
may  be  repeated  iteratively  until  convergence  is  met  within  a  prescribed  tol¬ 
erance. 

Criteria  for  selecting  an  appropriate  normalized  time  increment  (At/T) 
should  not  be  based  solely  upon  the  value  of  P  ,  i.e.,  the  implied  assumption 
of  acceleration  between  time  steps,  but  should  also  depend  upon  the  duration 
of  the  load.  Sufficiently  small  time  steps  must  be  used  to  provide  an  approp¬ 
riate  description  of  the  loading  function. 

As  an  example  of  the  numerical  parameters  which  might  be  used:  for  lin¬ 
ear  acceleration,  p  =  1/6  ;  the  normalized  time  increment  limit  for  convergence 
=  0.389;  and  the  normalized  time  increment  limit  for  stability  =  0.551  (Refer¬ 
ence  4).  Also,  the  time  increment  typically  should  not  exceed  about  1/10  of 
the  duration  of  the  load. 
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Table  2.1.  Speicher-Brode  relationships  for 
1-MT  surface  burst  weapon. 


SO  -  Surface  overpressure  (psi) 

-  Range  (kft) 

-  Time  of  arrival  (ms) 

^  -  Positive  phase  duration  (ms) 

-  Positive  phase  impulse  (psi-ms) 


PSO 

R 

t 

a 

fcd 

I 

P 

1.0 

46.767 

34,939.3 

4,376.8 

1,955 

2.0 

28.025 

18,646.6 

3,930.7 

3,285 

3.0 

21.122 

12,797.8 

3,583.2 

4,301 

4.0 

17.505 

9,817.7 

3,308.3 

5,136 

5.0 

15.248 

8,009.4 

3,085.3 

5,852 

6.0 

13.685 

6,790.6 

2,900.1 

6,482 

7.0 

12.527 

5,910.1 

2,743.4 

7,043 

8.0 

11.627 

5,242.1 

2,608.7 

7,548 

9.0 

10.903 

4,717.0 

2,491.5 

8,004 

10.0 

10.304 

4,292.4 

2,388.5 

8,415 

20.0 

7.264 

2,310.1 

1,785.3 

10,596 

30.0 

6.009 

1,609.0 

1,519.1 

11,366 

40.0 

5.282 

1,245.9 

1,378.0 

12,741 

50.0 

4.794 

1,022.5 

1,297.1 

14,598 

60.0 

4.437 

870.6 

1,249.3 

16,469 

70.0 

4.161 

760.2 

1,221.5 

17,986 

80.0 

3.939 

676.2 

1,206.2 

18,969 

90.0 

3.755 

610.0 

1,199.1 

19,482 

100.0 

3.600 

556.4 

1,197.5 

19,743 

200.0 

2.752 

305.0 

1,274.3 

28,198 

300.0 

2.366 

215.1 

1,357.8 

35,047 

400.0 

2.131 

167.9 

1,418.5 

38,823 

500.0 

1.967 

138.6 

1,461.9 

42,076 

600.0 

1.844 

118.4 

1,493.2 

45,206 

700.0 

1.746 

103.7 

1,516.2 

48,225 

800.0 

1.666 

92.4 

1,533.3 

51,112 

900.0 

1.599 

83.4 

1,546.0 

53,856 

1,000.0 

1.541 

76.1 

1,555.5 

56,460 

2,000.0 

1.214 

41.5 

1,575.7 

76,529 

3,000.0 

1.058 

28.9 

1,557.8 

90,071 

4,000.0 

0.960 

22.3 

1,535.2 

100,208 

5,000.0 

0.891 

18.3 

1,513.2 

108,274 

6,000.0 

0.838 

15.5 

1,493.1 

114,954 

7,000.0 

0.796 

13.4 

1,474.9 

120,646 

8,000.0 

0.761 

11.9 

1,458.1* 

125,602 

9,000.0 

0.732 

10.7 

1,443.4 

129,994 

10,000.0 

0.707 

9.7 

1,429.7 

133,940 

i ... — i 
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FORCING  FUNCTION 


RESISTANCE  FUNCTION 


Figure  2.1.  SDOF  analytical  model 
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a.  Nuclear  overpressure  decay  curve. 


b.  Normalized  pressure  history. 

2.3.  Typical  overpressure  decay  curves  for  a  specific 
nuclear  weapon  at  various  overpressures. 
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CHAPTER  3 


SOLUTIONS 


3.1  RESPONSE  CHARTS 

Once  the  differential  equations  have  been  nondimensionalized,  the  solu¬ 
tions  can  be  plotted  over  specific  ranges  of  the  nondimensionalized  parameters 
to  produce  a  practical  design  tool.  Normalized  response  charts  for  triangular 
loads  (Reference  2)  and  bilinear  loads  (Reference  5)  are  available  in  the  open 
literature. 

Charted  solutions  for  the  Speicher-Brode  overpressure  description  are 
provided  in  Figures  3. 1-3.9.  It  should  be  repeated  that  because  the  normal¬ 
ized  pressure  history  definition  is  a  function  of  two  independent  variables, 
i.e.,  time  and  peak  overpressure,  a  unique  solution  to  the  differential  equa¬ 
tion  cannot  be  obtained.  However,  uniqueness  was  obtained  in  this  case  by 
holding  the  peak  overpressure  constant,  thereby  creating  an  equation  with  a 
single  independent  variable.  As  a  result,  separate  charts  were  plotted  for 
each  peak  overpressure. 

3.2  SOLUTION  ACCURACY 

Response  charts  are  provided  for  peak  overpressures  ranging  from  10  psi 
to  50,000  psi.  For  design  purposes,  reasonable  accuracy  may  be  obtained  at 
intermediate  overpressures  by  interpolating  between  the  given  charts.  If  more 
accurate  solutions  are  desired,  a  computer  code  is  provided  in  Appendix  A  that 
will  give  numerically  accurate  solutions  for  any  desired  overpressure. 

A  comparison  of  solutions  obtained  from  the  response  charts  with  solu¬ 
tions  obtained  from  the  computer  code  has  indicated  that  logarithmic  interpo¬ 
lation  formulas  will  in  general  provide  more  accurate  solutions  between  re¬ 
spective  charts.  Although  linear  interpolation  between  any  of  the  charts  will 
provide  reasonable  answers,  more  accurate  solutions  can  sometimes  be  obtained 
with  logarithmic  interpolation  formulas.  Appropriate  interpolation  formulas 
are  presented  in  Table  3.1. 

3.3  EXAMPLE  PROBLEMS 

The  number  of  ways  that  the  response  charts  may  be  used  are  too  numerous 
to  sufficiently  describe  here.  In  addition  to  the  fact  that  the  user  may 
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enter  the  charts  with  several  different  coabinations  of  input  parameters , 
there  are  aany  different  ways  to  arrive  at  those  parameters  themselves,  e.g., 
natural  period,  yield  deflection,  and  maximum  resistance.  Hence,  the  follow¬ 
ing  examples  are  only  provided  to  illustrate  the  actual  mechanics  of  using  the 
charts.  It  is  assumed  that  the  user  will  know  how  to  compute  the  basic  prop¬ 
erties  of  the  physical  model. 

In  this  example,  the  maximum  response  is  determined  for  a  system  with  a 
natural  period  of  0.300  seconds,  a  maximum  resistance  of  10  psi,  and  a  peak 
elastic  deflection  of  0.5  inch.  The  structure  is  subjected  to  a  20-KT  surface 
burst  weapon  at  25  psi. 

In  the  response  curves,  the  natural  period  of  the  structure  has  been  nor¬ 
malized  with  respect  to  the  positive  phase  duration  of  the  load.  To  find  the 
positive  phase  duration  for  the  prescribed  loading,  see  Figure  2.2  (t^  =  1.6 
seconds) . 

By  cube-root  scaling,  the  duration  of  a  20-KT  weapon  is  coaqmted  as 


"d 

1.6  s 


20  KT  ' 
,000  KT j 


t^  =  0.434  s 


Computing  the  response  chart  input  parameters 


0.434 

0.300 


=  1.45 


—  =  =  2.50 

r  10 
m 


Now,  entering  the  10-psi  response  chart  (Figure  3.1) 


Mio  =  r  =  16 


and,  entering  the  50-psi  response  chart  (Figure  3.2) 


H50  =  6.5 
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Interpolating  (Table  3.1) 


M0c  =  exp 


In  /25> 

ln(50, 


(8) 


■  10  (rs) 


+  ln(6.5)  =  9.58 


The  predicted  maximum  response  is 


x  =  U-.x  =  9.58  (0.5  inch)  =4.79  inches 
m  r25  y 


Similarly,  the  estimated  time  to  maximum  response  is 


t  =0.38  second 
m 


The  maximum  response  as  determined  directly  by  computer  analysis  is  indi¬ 
cated  in  Appendix  B.  In  this  case,  there  is  an  interpolation  error  of  approx¬ 
imately  2  percent.  Using  a  linear  interpolation,  the  maximum  response  would 
be  predicted  as  p  =  12.4  ,  which  is  an  error  of  approximately  32  percent. 

In  a  second  example,  consider  a  system  with  a  natural  period  of  0.50  sec¬ 
ond  and  a  maximum  resistance  of  50  psi.  Using  a  75-KT  weapon,  the  peak  over¬ 
pressure  for  which  the  structure  can  maintain  a  ductility  of  15  may  be  ob¬ 
tained  by  the  following  procedure. 

Assume  a  peak  overpressure  of  100  psi.  From  Table  2.1,  the  positive 
phase  duration  for  a  1-MT  weapon  is  1.198  seconds.  For  a  75-KT  weapon,  the 
duration  is  scaled  as  follows: 


vd 

1.198  s 


(  75  KT  \ 
■  \l,000  KT/ 


td  =  0.505  s 

Now,  computing  the  nondimensionalized  parameters  based  on  the  assumed 
overpressure 

=  ^505  =  j 

T  0.050  10,1 


—  =  ^  =2.0 
r  50 
m 
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Using  the  100-psi  response  chart,  a  ductility  of  28  is  read.  That  is 
higher  than  our  Uniting  ductility.  Therefore,  try  a  lower  peak  overpressure 
of  50  psi  (Figure  3.2).  The  positive  phase  duration  of  1.297  seconds  for  a 
1-MT  weapon  at  50  psi  scales  to  0.547  second  for  a  75-KT  weapon. 

The  new  input  paraneters  are 


0.547 

0.050 


=  10.9 


=  §0  =i.o 
50 


Fron  the  50-psi  response  chart  (Figure  3.2),  a  much  lower  ductility  of 
3.2  is  read.  Finally,  interpolating  between  the  50-  and  100-psi  values  (Table 
3.1) 

ps  *  jzh  xtj-  •  <100  -  50>  +  50 

P  =  74  psi 
s  r 

The  correct  solution  was  determined  by  computer  analysis  to  be  approxi¬ 
mately  77  psi.  The  interpolation  error  is  approximately  4  percent. 

3.4  COMPUTER  CODE 

The  code  provided  in  Appendix  A  was  written  in  Fortran  for  the  Tektronix 
4081.  It  is  an  interactive  code  that  has  been  modularized  to  minimize  memory 
requirements.  The  code  should  be  easily  adaptable  to  other  machines  with  mod¬ 
ification  or  elimination  of  the  plot  routines  and  overlay  structure.  The  user 
should  also  find  it  quite  simple  to  add  alternate  loading  definitions  to  the 
code,  if  so  desired. 

Upon  initialization  of  the  program,  a  list  of  the  available  user  options 
is  displayed  to  the  screen.  These  options  and  a  brief  description  of  each 
follow.  A  sample  session  is  included  in  Appendix  8. 

X  -  EXIT  PROGRAM  -  The  user  may  exit  the  program  any  time  there  is  a 
request  to  SELECT  OPTION. 


iw. 
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ANALYZE  -  The  Newmark  0  analysis  is  performed  and  results  are  dis¬ 
played.  A  complete  response  history  including  normalized  displace¬ 
ments,  velocities,  and  accelerations  may  be  optionally  displayed. 
DISPLAY  DATA  -  All  numerical  parameters,  the  current  loading  des¬ 
cription,  and  structural  parameters  are  displayed  to  the  screen. 
HELP  -  A  list  of  the  available  options  is  provided  upon  request. 
INITIAL  PASS  -  This  option  assists  the  new  user  in  making  a  first 
run  through  the  program.  The  experienced  user  may  find  it  to  be  a 
beneficial  option  when  no  modifications  in  the  numerical  parameters 
are  required  prior  to  analysis. 

LOADING  FUNCTION  -  The  user  may  select  a  Speicher-Brode  or  triangu¬ 
lar  loading  description.  The  weapon,  peak  overpressure,  and  height 
of  burst  are  requested  if  the  Speicher-Brode  description  is 
selected.  Peak  overpressure  and  duration  are  requested  for  the 
triangular  function. 

NUMERICAL  PARAMETERS  -  Certain  numerical  parameters  may  be  changed 
by  the  user.  The  values  of  Beta,  the  integration  increment,  the 
convergence  tolerance,  and  the  number  of  pressure  stations  to  be 
computed  may  be  changed  by  selecting  this  option.  Default  values 
are  displayed. 

OUTPUT  DEVICE  -  When  large  quantities  of  data  are  to  be  output,  the 
user  may  select  for  that  data  to  be  routed  to  the  printer. 

PLOT  -  Normalized  displacements,  velocities,  accelerations,  and  a 
description  of  the  loading  function  may  be  displayed  in  graphical 
form  upon  request. 

SET  VALUES  FOR  ANALYSIS  -  The  natural  period  of  the  structure  and 
its  peak  resistance  are  requested  by  this  option. 

TITLE  -  Output  data  and  plots  display  a  title  that  may  be  changed 
by  selecting  this  option. 
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Table  3.1.  Interpolation  formulas. 


Type  of 
Interpolation 


Linear 


Natural 

logarithmic 


Notation:  P 

s 


To  ' 

Find  _  Formula 


Suggested 
Ranges  (psi) 


Mp  .  (P»  ’  P.l> 


50-100 


(Mp  Mi) 

P-  P-  =  tP-  - 4-  *  (P.0  -  Po1)  +  P, 


s  '  (M2  -  Mj) 


s2  si'  si 


M 


P 

s 


10-50 

100-200 

200-500 

500-1,000 

1,000-10,000 

10,000-20,000 

20,000-50,000 


-  overpressure  for  which  solution  is  desired 

-  overpressure  of  response  chart  before  Pg 

-  overpressure  of  response  chart  after  Pg 

-  ductility  at  desired  overpressure 
ductility  from  response  chart  before  Pg 
ductility  from  response  chart  after  Pg 


Figure  3.1.  Maximum  response  of  an  elastic-plastic  system 
subjected  to  a  10-psi  Speicher-Brode  pulse. 
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Figure  3.6.  Maximum  response  of  an  elastic-plastic  system 
subjected  to  a  1,000-psi  Speicher-Brode  pulse. 
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Figure  3.7.  Maximum  response  of  an  elastic-plastic  system 
subjected  to  a  10,000-psi  Speicher-Brode  pulse 
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CHAPTER  4 


CONCLUSIONS 

The  analysis  and  design  of  structures  to  resist  loadings  from  nuclear 
weapons  involve  many  parameters  that  are  difficult  to  define.  However,  if 
those  unknowns  can  be  determined  or  at  least  approximated,  then  a  straightfor¬ 
ward  solution  can  be  readily  obtained. 

Two  different  ways  for  arriving  at  that  solution  for  elastic-plastic  sys¬ 
tems  with  simulated  nuclear  loading  definitions  have  been  provided  in  this 
paper.  The  computer  solution  offers  a  quick  and  accurate  way  to  investigate 
the  effects  of  different  parameters  on  the  response.  However,  if  a  computer 
is  not  readily  available  or  if  a  rough  approximation  to  the  numerically  accu¬ 
rate  solution  is  adequate,  then  the  response  charts  are  most  beneficial.  The 
response  charts  provide  an  efficient  method  for  preliminary  design  calcula¬ 
tions  and  for  parametric  investigations  into  such  things  as  the  effect  of  peak 
overpressure,  weapon  yield,  structural  resistance,  and  structural  stiffness  on 
the  maximum  response  of  a  structure. 

The  charts  do  provide  numerically  accurate  solutions  for  the  specified 
overpressures,  but  interpolation  to  those  overpressures  not  included  on  a 
chart  does  introduce  some  error  into  the  analysis. 
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APPENDIX  A 

COMPUTER  CODE  LISTING 
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C******************************************************************** 


C  PROGRAM  FOR  DETERMINING  THE  NORMALIZED  RESPONSE  OF  AN  SDOF  * 

C  ELASTIC-PLASTIC  MODEL  SUBJECTED  TO  A  SPE I CHER- ERODE  OR  * 

C  A  TRIANGULAR  PULSE.  RESPONSES  ARE  NORMALIZED  BY  THE  ELASTIC  * 

C  DEFLECTION.  LES  GUICE  AUGUST  1983  * 

C******************************************* ************************* 
IMPLICIT  INTEGER*2< I-N, *  > 

INTEGER*1  TITLE<72) 

DIMENSION  P(250) ,RNU<501> »RNUD<501> ,RNUDD<501> 


DATA  TITLE/'  ','N', 

1  '  R' »'E' » 'S' » 'P' » 

2  '  ' » 'S'» 'D' » 'O' » 'F'» 

3  'P' » 'L' » 'A't 'S' »'T' » 

4  '  ' »  '  ' » '  '  »  '  '  , '  ' » 


'O' » 'R' » 'M' » 
'O' » 'N' » 'S' , 
'  ' ,  '  E '  »  '  L '  , 
'I't'C'f'  ', 


'A'  , 
'E'  » 
'A'  , 
'M'  , 


'L't'I'.'Z'f 'E't 
'  ' » 'O' » 'F' » '  '» 
'S' . 'T' » 'I' » 'C' » 
'O'.'D'f 'E'»'L'f 


'D' 

'A' 
/  ^  / 

/  / 
/  / 


NOPT=ll 

C 

C  INITIALIZE  ALL  VALUES 
C 

CALL  GRSTRT (9) 


CALL  BFACT  <0, 'USRJANLOLY .LIB  ') 


1  LU0=4 

J0F=4 
IPASS=0 
NSTP=501 
NST=500 
NT-250 


NP=250 
W=1000.0 
FKT= 100.0 


HOBK-O • 0 
BT=0. 1&&67 
T0L=1 . 0E-3 


TI 1=0 . 01 
TI2=0.01 
TN=1 .0 


RMX=1 .0 

10  CALL  OVLINM 'HELPOP  '»LUO> 

C 

98  CALL  OVLINK< 'GETOFT  ',NOF'T,JOF> 

GO  T0< 99999, 1000, 2000 t 1 0 , 3000 , 3000 , 4000, 5000 , 6000, 
1  7000,8000) ,JOP 


C 

C  ANALYZE  BY  NEUMARK  BETA  METHOD 
C 

1000  CALL  OVLINK< 'NEWMRK  ' ,LUO*P,NT,NP,RNU,RNUD,F:NUDD,NSTP,NST» 

1  TD , TN , T I 1 , T I 2 , PKT , RMX , BT , TOL , NNU , T I NC ) 

IPASS=2 
GO  TO  98 
C 

C  DISPLAY  DATA 
C 

2000  CALL  OVLINK< 'DISPOP  ' » LUO , JOF , LF UNC , W, PKT , TD , T A , R» HOBK , TN,RMX , 
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1  DPT  »P»NT.NP»PKT.  TITLE  * BT f TI 1 , TI2. TOL » NP > 

IF<J0P.EQ.5)60  TO  7000 
GO  TO  98 
C 

C  LOADING  FUNCTION 
C 

3000  CALL  OVLINM 'LOADOP  ' » LUO»  JOP»LFUNC » W»  PKT »TD»  TA»  RfHOBK » DPT » 
1  P» NT  f  NP) 

IPASS=1 

IF<J0P.EQ.5)G0  TO  2000 
GO  TO  98 

C 

C  NUMERICAL  PARAMETERS 

C 

4000  CALL  OVLINM '  NUMOF'  ' t LUO » BT  t  T 1 1 »  T 1 2 » TOL » NP ) 

GO  TO  98 

C 

C  OUTPUT  DEVICE 

C 

5000  CALL  OVLINM 'OUTOP  ' »LUO) 

GO  TO  98 

C 

C  PLOT  RESULTS 
C 

6000  CALL  OVLINM 'PLOTOP  '  » LUO »  T I TLE  t RNU» RNUD » RNUDD » NSTP » 

1  Pr NT  t NP »  TD »  NNU  *  TINC» I PASS ) 

60  TO  98 

C 

C  SET  VALUES  FOR  ANALYSIS 

C 

7000  CALL  OVLINM 'SETOP  ' . LUO . TN » RMX ) 

IF  <  JOP .EG • 5)60  TO  1000 
GO  TO  98 
C 

C  TITLE 
C 

8000  CALL  OVLINM 'TITLOP  ' »  TITLE ) 

60  TO  98 

C 

C  EXIT  PROGRAM 

C 

99999  STOP 
END 


SUBROUTINE  HELPOP(LUO) 

IMPLICIT  INTEGER4t2( I-N,#) 

CALL  VTSCA(3) 

CALL  HIBRN8( 1 ) 

WRITE(LUO»4t>  '  ' 

WRITE (LUO,  40 '  NORMALIZED  RESPONSE  OF  AN  S  DOF  ELASTIC-'. 

1  'PLASTIC  MODEL' 

URITE<LUO»*>'  ' 

WRITE ( LUO. 4c ) '  NOTE {  DEFAULT  UNITS  ARE  INCHES.  POUNDS.  AND  SEC' 
WRITE  (LUO.  4c) '  UNLESS  OTHERWISE  SPECIFIED.' 

WRITE  (LUO.  4c) '  ' 

WRITE  (LUO.  4( ) '  AVAILABLE  OPTIONS: ' 

WRITE(LUO»4c) '  X  -  EXIT  PROGRAM' 

UR  I TE  (  LUO  .  4c )  '  A  -  ANALYZE  ' 

WRITE(LU0.4c)  '  D  -  DISPLAY  DATA' 

WRITE(LU0.4C)  '  H  -  HELP' 

WRITE(LU0.4c) '  I  -  INITIAL  PASS' 

WRITE  (LUO.  4c) '  L  -  LOADING  FUNCTION' 

WRITE ( LUO. 4c )  '  N  -  NUMERICAL  PARAMETERS' 

WRITE(LU0.4c)  '  0  -  OUTPUT  DEVICE' 

WRITE ( LUO . 4c )  '  P  -  PLOT' 

WRITE ( LUO. 4t )  '  S  -  SET  VALUES  FOR  ANALYSIS' 

WRITE (LUO .40  '  T  -  TITLE' 

WRITE ( LUO.  4c ) '  ' 

RETURN 

END 


SUBROUTINE  GETOPT ( NOPT . JOP ) 

IMPLICIT  INTEGER4C2(I-N.*) 

INTEGER4C2  JOPT(ll) 

DATA  JOPT/'X  ' . ' A  '»'D  '»'H  '.'I  '»'L  '.'N  '.'0  '»'P  ', 

1  'S  ' .  'T  '/ 

WRITE  (  6 . 4c )  '  SELECT  OPTION!  <  X.  A.D.H.  I  .L.N.O.P.S.T)  ' 

CALL  READC ( I OPT .0.0) 

WRITE  (6. 4c)  '  ' 

DO  100  1=1, NOPT 
IFdOPT.EQ.  JOPT(I)  >JOP=I 
RETURN 
END 


o  non 


SUBROUT I NE  NEWMRK ( LUO » P *  NT  * NP  *  RNU » RNUD , RNUDD » NSTP  *NST * 

1  TD*TN*TI1 *TI2*PKT  * RMX*BT  *TQL*NNU*TINC> 

THIS  SUBROUTINE  IS  USED  TO  INTEGRATE  THE  SDOF  EQUATIONS 
OF  MOTION  BY  THE  NEWMARK-BETA  METHOD 

IMPLICIT  INTEGER*2<I-Nf*> 

DIMENSION  P < NT ) * RNU < NSTP ) * RNUD ( NSTP  >  t RNUDD < NSTP ) 

DATA  INO/'N  '/ 

PI=3. 1415927 
ANUM=4.0*PI**2.0 
C  COMPUTE  MAXIMUM  RESPONSE  FOR  THE  RATIO  OF  TD/TN  AND  PKT/RMX 
POR=PKT/RMX 

C  DETERMINE  TIME  INCREMENT  TO  SATISFY  CONVERGENCE  AND 
C  PROPER  DESCRIPTION  OF  LOADING  FUNCTION 
TI=TI 1 
TINC=TI*TN 

IF (TD.GT ♦TN)GO  TO  1000 

TI=TI2 

TINC=TI*TD 

C  SET  INITIAL  VALUES  OF  NORMALIZED  DISPLACEMENT (RNU) » 

C  VELOCITY(RNUD) t  AND  COMPUTE  INITIAL  ACCELERATION ( RNUDD ) 

1000  RNU< 1 )=0.0 

RNUMX=RNU( 1 ) 

TMXT=0.0 
RNUD( 1 ) =0 . 0 
RNUDD  <  1  >  ~ANIJM*P0R 
TIDTN=TINC/TN 

C  PERFORM  THE  INTEGRATION  ITERATIVELY 

C  BY  ASSUMING  AN  ACCELERATION  AND  CHECKING  CONVERGENCE 
C 

DO  1500  J=1 »NST 
JJ=J+1 

RNUDD  <  J J ) =RNUDD ( J ) 

C  DETERMINE  THE  NORMALIZED  PRESSURE  THIS  TIME  STEP 
CALL  PRESS ( P  *  NT  » NP » TD » J J »  T INC » SI » FT ) 

C  COMPUTE  NORMALIZED  DISPLACEMENT  AND  VELOCITY 

1100  RNU< JJ)=RNU< J)+TIOTN*RNUD( J ) +TIDTN*TIDTN* ( 0 . 5-BT > *RNUDD( J)+ 
1  TIDTN*TIDTN*BT*RNUDD< JJ) 

RNUD( JJ ) =RNUD (J)+0.5*TIDTN#( RNUDD <  J ) +RNUDD ( JJ ) ) 

TNU=RNU ( JJ  > 

C  DETERMINE  IF  RESPONSE  IS  IN  PLASTIC  RANGE 
IF(RNU( JJ) . GT . 1 4  0 ) TNU-l . 0 
C  CHECK  CONVERGENCE 

RNUCN=ANUM*(FT*POR-TNU> 

IF(ABS(RNUCK-RNUDD(  JJ)  )  .LE.TODGO  TO  1200 
IF (RNUCK.EQ. RNUDD ( J) >G0  TO  1200 
C  NOT  CONVERGED ,  ASSUME  NEW  ACCELERATION 
RNUDD <JJ)=RNUCK 
GO  TO  1100 

C  CONVERGED*  DETERMINE  IF  A  MAXIMUM 


1200  IF<RNU< JJ) . LE.RNU< J> >G0  TO  2000 
RNUMX=RNU< JJ) 

TM=SI 

FM=FT 

NNU=JJ 

1500  CONTINUE 

URITE<LUOf #) '  caution:  SOLUTION  NOT  CONVERGED ! ' 

c 

c  DISPLACEMENT  IS  DECREASINGf  MAX.  AT  PREVIOUS  STEP 
C  URITE  OUT  RESULTS  OF  ANALYSIS 
2000  WRITE(LUO»#) ‘  ' 

WRITE (LUOf*) '  MAXIMUM  DUCTILITY  ='»RNUMX 

URITE ( LUO » # ) '  OCCURS  AT  TIME  ='»TM 

WRITE < LUO » # ) '  NORMALIZED  PRESSURE  AT  THIS  TIME 

WRITE ( LUO  t%) ‘  ' 

WRITE(LUOrft) '  NATURAL  PERIOD  =  ',TN 

WRITE<LUO»#) 7  LOAD  DURATION  =  ',TD 

WRITE(LUO»#) '  TIME  INCREMENT  =',TINC 

WRITE ( LUO » * ) '  ' 

WRITE < 6 » # > '  COMPLETE  HISTORY?  CY/N):' 

CALL  RE ADC ( I ANS  >  0 » 0 ) 

WRITE(6»*) 7  ' 

IF <IANS.EQ.INO> RETURN 
WRITE(LUO»*) '  ' 

WRITE ( LUO >  * )  '  NOTE.*  RESPONSES  ARE  NORMALIZED." 

WRITE(LUO>*) '  ' 

WRITE<LUO.  * )  •'  TIME  DISP  VEL 

T=0.0-TINC 

DO  2500  1=1 »JJ 
T=T+TINC 

2500  WRITE* LUO , 251 0  >  T , RNU  < I > • RNUD  < I ) » RNUDD ( I ) 

2510  F  ORMAT  <  2  <  2X » F 1 0 . 5 ) «2(2X»F10.5> ) 

WRITER  LUO »  #  )  '' 

RETURN 

END 


cnn  n  nnnnnnn  non 


SUBROUT I NE  PRESS  <  P , NT , NP , TD » J J »  T I NC , SI , FT ) 

SUB  FOR  COMPUTING  NORMALIZED  PRESSURE  AT  A  GIVEN  TIME  STEP 

IMPLICIT  INTEGER*2<I~N»*> 

DIMENSION  p < NT ) 

SI  =  REAL  TIME  FOR  THIS  STEP 
SID  =  NORMALIZED  TIME  AT  END  OF  DURATION 
FT  =  NORMALIZED  PRESSURE  AT  GIVEN  TIME 
TINC  =  TIME  INCREMENT 

DELTT  =  TIME  INCREMENT  BASED  ON  NUMBER  OF  PRESSURE  STATIONS 
SI=TINC*FLOAT< JJ-1 ) 

IF  TIME  IS  BEYOND  POSITIVE  DURATION,  SET  PRESSURE  =  0 
IF  <  SI . LT . TD ) GO  TO  500 
FT =0 • 0 
RETURN 

FIND  PRESSURE  STATION  ON  EITHER  SIDE  OF  GIVEN  TIME 
DELTT  =  TD/FLOAT<NP-:l  ) 

1ST 1=INT2( SI /DELTT >+I 
ST 1=FL0AT ( IST1 ) 

ST2=ST1  +  1 .0 
IST2=IST1+1 
P11=P< IST1 ) 

PI 2=P  < IST2 ) 

TST1=(ST1-1 .0 >*DELTT 
TST2=TST 1+DELTT 
C  INTERPOLATE  TO  FIND  PRESSURE  THIS  TIME  STEP 
FT  =  (  PI  1--P1 2  >#<SI -1  ST2  )  /  <  T ST  l  -  rST2  )  +P1 2 
RETURN 
END 
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C 

c 

c 


2010 


SUBROUT I NE  D I SPOP  ( LUO » JOP » LFUNC » U » PKT » TD  »  T A  » R » HOBK , TN » RMX  . 
1  DPT .  F'.NT  .NF’.PKT  .TITLE.BT  .TI1.TI2.T0L.NP) 

IMPLICIT  INTEGER*2(I-N»#> 

INTEQER*2  IFU(3> 

INTEGERS  TITLE(72> 

DIMENSION  P<NT> 

DATA  INO/'N  V.IFU/'S  '.'T  '.'X  '/ 

DISPLAY  DATA 

WRITE(LUO»*) '  ' 

WRITE (LUO. 2010 ) (TITLE (I).I=1» 72) 

FORMAT  < IX*  72 A 1 ) 


2300 

2310 


=  ' »TI1> 
=  ' » TI2 1 


2050 

2100 


WRITE(LUO. *> '  ' 

WRITE ( LUO.  # ) ' 

WRITE  <  LUO  » # ) ' 

WRITE < LUO t 4) ' 

WRIT£(LUO»*> ' 

WRITE(LUO.*> ' 

WRITE(LUO.*) ' 

WRITE(LUO»*) '  ' 

WRITE(LUO»*) ' 

IF  <  LFUNC . EQ . IFU< 1 ) ) QO 
IF  < LFUNC. EQ. IFU(2) )G0 
GO  TO  2200 

WRITE(LUQ»*) '  WEAPON  (KT)='»W 

WRITE (LUO.*) '  PEAK  OVERPRESSURE  =  ',F'KT 

WRITE (LUO .* ) '  POSITIVE  PHASE  DURATION  =  '.TD 

IF (LFUNC. EG. IFU<2) >G0  TO  2200 


NUMERICAL  PARAMETERS  J ' 

BETA  =  ' .  BT 
INTEGRATION  INTERVAL 
INTEGRATION  INTERVAL 
TOLERANCE  ='.TOL 
NUMBER  OF  PRESSURE  STATIONS  =' 

LOAD  description:' 

TO  2050 
TO  2100 


OF 

OF 

.NP 


PERIOD' 

DURATION' 


TIME  OF  ARRIVAL  =  '»TA 
RANGE  <  KFEET  >  =  '.R 
HEIGHT  OF  BURST  (KFEET)  =' 
2250 


2200 


STRUCTURE  DESCRIPTION! 
NATURAL  PERIOD  =  ',TN 


FREQUENCY  * ' .FREQ 
PEAK  RESISTANCE  =' 


RMX 


2250 


WRITE(LUO»*) ' 

WRITE ( LUO.  * ) ' 

WRITE(LUQ»* ) ' 

I F  (  JOF'  •  EG  .  5 )  GO  TO 
WRITE(LUO.*> '  ' 

WRITE(LUO,*) ' 

WRITE(LUO»*) ' 

FREQ=1 . O/TN 
WRITE(LUO.*> ' 

WRITE (LUO.*) ' 

WRITE(LUO»*) '  ' 

WRITE (6.*) '  PRESSURE-TIME  HISTORY?  (Y/N) 
CALL  RE ADC (IANS. 0.0) 

WRITE(6»*) '  ' 

IF( IANS. EG. I NO) RETURN 
T=0.  O-DF'T 
WRITE(LUO»*> '  ' 

WRITE ( LUO » # ) '  TIME  PRESSURE 

WRITE(LUO.*> '  ' 

DO  2300  I=1,NP 
PTPKT=P(I)*PKT 
T=T+DFT 

WRITE(LU0.2310)T .PTF’KT .P(  I ) 

FORMAT ( 3 ( 2X . F 1 0 . 4  )  ) 

RETURN 

END 


.HOBK 


NORMALIZED' 
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SUBROUT 1 NE  LOADOP <  LUO . JOP » LFUNC » U  . PKT , TD » T A »  R . HOBK ,  DPT » 
1  P.NT.NP) 

IMPLICIT  INTE6ER*2<I-N.#> 

INTEGERS  IFU(3) 

DIMENSION  P ( NT  > 

DATA  IFU/'S  '  »  'T  '.'X  '/ 

LOADING  FUNCTION 

WRITE<6»*)'  SELECT  SPEICHER-BRODE .  TRIANGULAR. ' , 

1  '  OR  EXITJ  (S.T.X)' 

CALL  READC( LFUNC. O.O) 

WRITER.*)'  ' 

DO  3100  1=1.3 
IF (LFUNC . EQ . IFU< I ) ) IFUOP=I 
3100  CONTINUE 

GO  T0<3200. 3300. 3500) . 1FU0P 
3200  URITE<6,#)'  ENTER  WEAPON(KT>:/ 

READ(5.*)W 

WRITE ( 6 » # > '  ENTER  PEAK  OVERPRESSURE:' 

READ (5.*) PKT 

WRITE<6.*)'  ENTER  HEIGHT  OF  BURST  <  KFEET ).‘ ' 

READ (5.#) HOBK 

CALL  BRODEP < P . TD . DPT . W . PKT . NT » NP . HOBK . R » TA > 

TD=TD/1000.0 
TA=TA/1000.0 
DPT=DPT/ 1000.0 
GO  TO  3400 

3300  WRITE (6,#) '  ENTER  PEAK  OVERPRESSURE:' 

READ (5.*) PKT 

WRITE (6.#) '  ENTER  DURATION:' 

READ (5 » # ) TD 

CALL  BIGGS ( P , NT . NP . TD . DPT . PKT ) 

3400  REAEK5.3410) IDUMMY 

3410  FORMAT  <  A2 ) 

3500  RETURN 
END 
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C*TITLE  FILED  IN  BRODEP 

SUBROUTINE  BRODEP (PEST.  DP.  DT »  W.  P»  NPTS.NP.HOBK.R.TA) 
IMPLICIT  INTEGER  *2(I-N»#) 

DOUBLE  PRECISION  XP.YP.DPP 


DIMENSION  PEST(NPTS) 

TA  ARRIVAL  TIME  (MSEC) 

RANGE  IN  KFEET 

DP  POSITIVE  PHASE  DUR.  (MSEC) 

DT  INCREMENT  BETWEEN  PRESSURE  STATIONS 
IONE  =  0 

W13  =  W**(l./3.) 

YP  =  HOBK/W13 
DTT  =  DT/W13 

XP  =  ((1.58  *  2.  *  U  /  P)  **  <1.  /  3 . >  >/U13 
10  CALL  PT (XP. YP.DTT .PEST ( 1 ) .DPP. IONE.NPTS.NP.TA) 

IF  ( ABS ( PEST ( 1 > /P-1 ♦ )  .LE.  5.E  -  5)  GO  TO  20 
CALL  PT ( XP* 1 . 0000 1  *  YP»  DTT  >EST » DPP .IONE. NPTS . NP . TA ) 

XP  =  XP  +  ( PEST ( 1 )  -  P)  *  0.00001  *  XP  /  ( PEST ( 1 )  -  EST) 
. GO  TO  10 
20  R  =  XP*W13 
I0NE=1 

CALL  PT ( XP . YP , DTT . PEST . DPP . IONE . NPTS . NP. TA > 

DP  =  DPP*U13 
7A  =  TA#W13 

NOTE  -  REDEFINITION  OF  DT  BELOW  -  LKG  7/83 
DT=DTT*W13 


RETURN 

END 


C*TITLE  FILED  IN  FT 

SUBROUTINE  PT<X.Y»DTT.P»DP.IP.NPT.NP.TAS> 

C 

C  THIS  IS  AN  IMPLEMENTATION  OF  THE  SPEICHER  AND  BRODE 
C  OCT  1981  PIT.X.Y)  NUC.  HISTORY 

C  REF.  PSR  NOTE  385  WITH  MOD  FOR  P  ABOVE  lOOOO  PSI 
C  WALKER  NOV. 17 . 1982 
C  INPUT  PARAMETERS 
C  X  =  RANGE  <KFT/W**l/3> 

C  Y  =  HOB  <KRT/W**l/3> 

C  T  =  TIME  AFTER  TIME  OF  ARRIVAL  (MSEC) 

C  IP  =  0  FOR  PEAK  PRESSURE  ONLY 
C  DTT  =  TIME  STEP  SCALED  <MSEC/W**l/3> 

C 

C  OUTPUT  PARAMETERS 
C 

C  P  =  PRESSURE*TT(PSI) 

C  DP  *  POSITIVE  PHASE  DURATION  <MSEC> 

C  ***NOTE  T  IS  CHANGE  IN  THIS  ROUTINE. 

IMPLICIT  INTEGERS  (I-N.#) 

IMPLICIT  DOUBLE  PRECISION  (A-H.O-Z) 

REAL  P.  DTT.  TAS 
DIMENSION  P(NPT) 

C 

XLEAST=1 • E-4 
YLEAST=l.E-4 
ZLEAST=1 . E-4 
ZM0ST=100. 

IF  <  X . LT • XLEAST )  X=XLEAST 
IF(  Y ,  LT..  YLEAST  )  Y=YLEAST 
R=DSQRT<X*X+Y*Y> 

Z=Y/X 

IF  <  Z . LT . ZLEAST )  Z=ZLEAST 
IF(Z.GT.ZMOST)  Z=ZMOST 
CALL  PPEAK (X.Y.R.Z. DEL  TAP ) 

P  < 1 )  =  DELTAP 
IF< IP.EQ.O)  RETURN 

XM= 1 70. *Y/<1 .+337.* Y**.25>+.914*Y*#2. 5 
U= ( . 543-21 . 8*R+386 . *R*R+2383 . *R*R*R ) *R**8 

UD=2.99E-14-< 1 .91E-10*R*R>+<1 .032E-6*R**4>-(4.43E-6*R**6>+ 

8  ( 1 . 028+ ( 2 . 087*R )  +  ( 2 . 69*R*R  > ) *R**8 

U=U/UD 
TA=U 

IF ( X . LT . XM ) GO  TO  200 

W=  ( 1 . 086-34 . 605*R+486 , 3*R*R+2383  .  *R*R*R )  *R**8 

WD=3 . 0137E- 1 3-  ( 1.21 28E-9*R*R >+<4.1 28E-6*R**4 )  -  ( 1 . 1 1 6E-5*R**6 >  + 

8  <1 . 632+2. 629*R+2.69*R*R)*R**8 

W-W/WD 

TA=<U*XM/X)+<W*< 1 .-XM/X) ) 

200  CONTINUE 

S  =  1,-1.1E10*<Y**7)/<1.+1.1E10*Y**7>-<2.441E-8*Y*Y/ 

8  <1.+9.E10*Y**7>  )*<1  ,/<4.41E-ll+X**10> ) 

F  =  <.01477*<TA**.75)/fl.+.005836*TA>+7.402E~5*<TA**2.5)/ 

*  <1 .+1 .429E~8*TA**4.75)-.216)*S  +  .  7076-3 . 077E-5*TA*TA*TA/ 

8  <1.+4.367E-5*TA*TA*TA> 

G  =  10.+<77.58~64.99*<TA**.125>/(1.+.04348*DSQRT<TA>>)#S 
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H  =  2.753+.05601*TA/<l.+1.473E-9*TA**5>+(.01769*TA/ 
l  <l.+3.207E-10*TA**4.25>-.03209*<TA**1.25>/<l.+9.914E-8* 

t  TA**4>-1 .6>*S 

DP  »( <1640700. +24629. *TA+416. 15*TA*TA>/< 10880. +619, 76*TA+TA*TA> )* 
*  < . 4+ . 001204* ( TA** 1 . 5  >  /  <  1 .  + . 001 559*TA** 1 . 5  >  +  < . 0426+ . 5486* 

X  <TA**.25)/<1 .+.00357*TA**1 .5) )*S) 

C*** 

C  NOTE  MODIFICATION  BELOW  - 

C  COMPUTE  TIME  STEP  BY  DIVIDING  THE  POSITIVE  PHASE  DURATION 

C  BY  THE  NUMBER  OF  DESIRED  TIME  STEPS 

C*** 

DTT=DP/FLOAT  <NP-1 ) 

DO  500  I  =  1 t NP 
T  =  TA+(DTT*( 1-1 ) ) 

B  =  (F*  <  TA/T)**G+  < 1 ♦ -F ) * ( TA/T  >  **H>  *( 1 • -<  T-TA)/DP) 

IF  (X.LT.XM  .OR.  Y.GT.0.38)  GO  TO  1000 
XE  =  3.039*Y/( 1 .+6.7*Y) 

E  =  DABS< <X-XM)/(XE-XM> ) 

IF  (E.GT.50. >E=50. 

D  =  . 23+583000. *Y*Y/(26667. +1000000. *Y*Y>+.27*E+( .5-583000. *Y*Y/ 

X  (26667. +1000000. *Y*Y) )*E**5 

A  =  (D-l. )*(1.-(E**20)/<1.+E**20) > 

DT  =  474. 2*Y*(X-XM) **1.25 
IF(DT.LT.l.E-4)  DT=l.E-4 
GA  =  (T-TA)/DT 
IF (GA. GT . 400. )  GA=400. 

V  =  1 . + ( 3 . 28E1 1* ( Y**6 ) / < 1 . + 1 . 5E12*Y**6. 75 ) ) *(GA*GA*GA/ 

X  <  6 . 13+GA*GA*GA ) ) *  < 1 . / < 1 . +9 . 23*E*E  >  > 

C=  ( < 1 .04-240.9*<X**4)/( 1 .+231 ,7*X**4) >*(GA**7)/( ( 1 .+.923* 

X  ( GA ) **8 . 5 ) « < 1 +A ) ) )*(1 .-( (T-TA)/DP)**8)*2.3E13*(Y**9>/<  1 .+ 

X  2.3E13*Y**9> 

P<I>  =  DELTAP*( 1 .+A)*(B*V+C) 

GO  TO  495 
1000  CONTINUE 

P(I>  =  DELTAP*B 

C*** 

C  NOTE  MODIFICATION  BELOW  - 

C  NORMALIZE  ALL  PRESSURES  BY  DIVIDING  BY  PEAK  PRESSURE 

C*** 

495  P( I )=P( I ) /DEL TAP 

500  CONTINUE 
TAS  =  TA 
RETURN 
END 
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C*TITUE  FILED  IN  PPEAK 

SUBROUT INE  PPEAK <X»Y»R»Z» DEL TAP ) 

C 

C  SUB  FOR  PEAK  OVERPRESSURE  BY  SPEICHER  AND  BRODE t SEPT . 1 980 
C  MOD  BY  WALKER  FROM  MOD  BY  SPIECHER  AND  BRODE  MARCH  1981 
C  RANGE  TO  0.1  PS I 
C 

C  INPUT  PARAMETERS 
C  X*RANGE  <KFT/KT**l/3> 

C  Y=HOB  <KFT/KT**l/3> 

C  R=SQRT<X**2+Y**2>  ( KFT/KT**1 /3 ) 

C  Z=Y/X 
C  OUTPUT 

C  DELTAP=PEAK  PRESSURE  <PSI) 

C 

IMPLICIT  INT£GER#2  (I-N»*) 

IMPLICIT  DOUBLE  PRECISION  <A-H»0-Z> 

C 

A  =  1.22-<3.908*Z*Z)/(1.+810.2*Z**5> 

B  =  2.321+<<Z**18>/(l.+l.tl3*Z**18))*6.195-(.03831*Z**17>/ 

8  < 1 , + . 0241 5*Z** 17)+. 6692/ < 1 . +41 64 . *Z**8 ) 

C  =  4.153-(1.149*Z**1B)/<1.+1.641*Z**18)-1.1/<1.+2.771*Z**2.5> 

D  =  -4.166+<25.76*Z**1.75)/<1.+1 , 382*Z**18 >+8. 257*Z/ < 1 . +3. 219*Z > 
E  =  l.-(.004642*Z**18>/<l.+.003886*Z**18> 

F  =  .6096+<2.879*Z**9.25>/< 1 . +2 . 359*Z**14 . 5) -17 . 15*Z*Z/ 

8  <1.+71.66*Z*Z*Z> 

G  =  1.83+5.361*Z*Z/<1 .+.3139*Z**6> 

H  =-<64.67*Z**5+.2905)/<l .+441 .5*Z**5>-1 .389*Z/< 1 . +49 . 03*Z**5 > + 

8  <8.808*Z**1.5)/< 1 .+154.5*Z**3.5>+( .0014*R*R/< 1 . - . 158*R+ . 0486* 

8  R**1.5  +.00128*R*R) )*( 1 ./< 1 .+2.*Y> ) 

C 

DELTAP* 10 . 47/ ( R**A ) +B/ ( R**C ) +D*E/ ( 1 . +F*R**G ) +H 

RETURN 

END 


n  o 
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SUBROUTINE  BIGGS<P.NT,NP» TD.DTT.PKT) 

GENERATE  NONDIMENSIONALIZED  PRESSURES  FOR  TRIANGULAR  PULSE 
IMPLICIT  INTEGER*2<I-N.*> 

DIMENSION  P(NT) 

DTT=TD/FLOAT (NP-1 ) 

P( 1 )=PKT/PKT 

DO  500  1=1 .NT 
P<I)=P<1)-<I-1>*DTT/TD 
RETURN 
END 


SUBROUT I NE  NUMOP < LUO . BT . T 1 1 , T 1 2 » TOL . NP ) 

C 

C  NUMERICAL  PARAMETERS 
C 

IMPLICIT  INTEGER*2(I-Nt*) 

INTE6ER*2  NUOP(S) 

DATA  NUOP/'B  '.'I  '»'T  '»'N  '.'X  '/ 

N0P=5 
NUM0PT=5 
WRITE ( 6» * ) '  ' 

WRITE (6.#) '  SELECT  NUMERICAL  PARAMETER.*  (B»I»T »N»X>' 
WRITE<6»*>'  B  -  BETA' 

WRITE ( 6  >  * ) '  I  -  INTEGRATION  INCREMENT' 

WRITE(6»*)'  T  -  TOLERANCE' 

WRITE<6»*>'  N  -  NUMBER  OF  PRESSURE  STATIONS' 

WRITE(6»*>'  X  -  EXIT' 

WRITE(6.*>'  ' 

CALL  READC( I OPT »0»0) 

WRITE(6»*>'  ' 

NUM0PT=5 

DO  10  1  =  1  > NOP 

10  I F  <  I  OF'T  .  EQ .  NUOP  ( I )  )  NUMOF'T = I 

GO  T0<  1 00 » 200 »  300. 400.500 ) » NUMOF'T 
C 

100  WRITE ( 6 » #) '  ENTER  BETAS  (LINEAR  ACCELERATION  =  0.16667)' 

READ(5»*)BT 
GO  TO  498 

200  WRITE (6.#) '  ENTER  TIME  INCREMENT/PERIODS  (TI/TN  =  0.01)' 

READ (5»*)TI1 

WRITE (6  >  * ) '  ENTER  TIME  INCREMENT/LOAD  DURATIONS  ( TI/TD  =  0.01) 
READ ( 5  »  * ) TI2 
GO  TO  498 

300  WRITE(6.*)'  ENTER  CONVERGENCE  TOLERANCES  (TOL.  =  IE-3)' 

READ(5»*)T0L 
GO  TO  498 

400  WRITE ( 6  »  #  )  '  ENTER  NUMBER  OF  PRESSURE  STATIONS.*  (MAX.  =  250)' 

WRITE  <  6.  #  > '  CAUTIONS  MUST  COMPUTE  NEW  LOAD  FUNCTION!' 

READ(5.*)NF' 

498  RE  AIK  5.  499) I  DUMMY 

499  FORMAT (A2) 

500  RETURN 
END 
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SUBROUTINE  OUTOP<LUO> 

IMPLICIT  INTEGER*2<I-N>#> 

INTEGER*2  JDEV<2> 

DIMENSION  P(250>  »RNU( 1001 >  >RNUD< 1001 >  r RNUDD (1001) 
DATA  JDEV/'T  '*'P  '/ 

OUTPUT  DEVICE 

5000  URITE<6?*> '  SELECT  TERMINAL  OR  PRINTER:  (TrP)' 
CALL  READC <  I DE  V » 0  >  0 ) 

WRITE(6r*> '  ' 

IFdDEV.EQ.  JDEV<1>  )L.U0=6 
IF< I DEV . EQ ♦ JDEV ( 2 ) ) GO  TO  5100 
RETURN 
5100  LU0=8 

CALL  ASSIGN('LF:  ' »  LUO  *  IS  TAT  > 

RETURN 


o  no 


SUBROUT I NE  PLOTOP <  LUO  »  T I TLE  »  RNU » RNUD » RNUDD »  NSTP  r 
1  P»NT » NP > TD»NNU > T INC » I PASS) 

IMPLICIT  INTEGER*2< I-N»*> 

INTEGERtl  TITLE<72) 

INTEGER*2  JPLT(5> 

DIMENSION  P(NT>  r RNU (NSTP) r RNUD (NSTP) r RNUDD (NSTP) 
DATA  JPLT/'D  '  '.'A  ','L  't'X  '/ 

PLOT  RESULTS 

WRITE(6»*>'  SELECT  PLOTi  (D»V.A»L,X)' 

CALL  READC( I PL T >0*0) 

WRITE(6.*)'  ' 

JPL*5 

DO  6100  1=1 .4 
6100  IF(IPLT.EQ.JPLT(I> )JPL=I 
IF  <  JPL . EQ . 5 ) RETURN 
IF< IPASS.LT < 1 )G0  TO  6200 
IF (IPASS.LT. 2. AND. JPL. LE. 3) GO  TO  6200 
DT=TD/(FLOAT (NP-1 ) )*1000.0 
TINCP=TINC*1000.0 

CALL  PLOT  <  JPL  >  NSTP , NT » NP » RNU » RNUD , RNUDD »  P » T I TLE » 
1  NNU.TINCPf DT ) 

RETURN 

6200  WRITE(6»*)'  ERROR?  DATA  FOR  PLOT  NOT  AVAILABLE ! ' 

RETURN 
END 


SUBROUTINE  PLOT ( JPL .NSTP.NT  *NP. RNU » RNUD » RNUDD » P » T I TLE » 

1  NNU » T INC » BELT) 

C  SUBROUTINE  FOR  PLOTTING  RESPONSES 
IMPLICIT  INTEGER#2< I-N*#) 

INTEGERS  TITLE  <  72 ) 

DIMENSION  P(NT ) .RNU(NSTP) » RNUD(NSTP) . RNUDD (NSTP) 

DATA  K1 1/1 1/ »  K22/22/ .  K30/30/,  K3A/36/.  TZ/O./ 

C 

IF ( JPL«LT *4) CALL  TITLX< 'TIME  <MSEC)'»K11) 

CALL  NEXTPG 

C  WRITE  TITLE  TO  PLOT 
CALL  VTSCA<5> 

CALL  HIBRN8( 1 ) 

WRITE<6»90)TITLE 
90  FORMAT ( 10X.72A1 ) 

CALL  VTSCA<3> 

CALL  0IBRN8< 1 ) 

C 

GO  TO  <100.20) .300.400) . JPL 

100  CALL  TITLYK 'NORMALIZED  DISPLACEMENT  <D/XY)',K30> 

CALL  QPLOTD (TZ.TINC. RNU » NNU > 

CALL  NEXTPG 
RETURN 

200  CALL  TITLYK 'NORMALIZED  VELO  CITY  <V*T/XY>  '.K30) 

CALL  QPLOTD <TZ.T I NC.RNUD. NNU) 

CALL  NEXTPG 
RETURN 

300  CALL  TITLYK 'NORMALIZED  ACCELERATION  <A*T*T/XY>  '.K3&) 
CALL  QPLOTD <  TZ , T INC . RNUDD . NNU ) 

CALL  NEXTPG 
RETURN 

400  CALL  TITLXt 'TIME  (MSEC) '.Kll) 

CALL  TITLYK 'NORMALIZED  PRESSURE  '»K2 2) 

CALL  QPLOTD (TZ.DELT.P.NP) 

CALL  NEXTPG 

RETURN 

END 


non 


SUBROUTINE  SETOP(LUO*  TN.RMX) 
IMPLICIT  INTEGER*2(I-N.#) 

SET  VALUES  FOR  ANALYSIS 

WRIT£<6»*)'  ENTER  NATURAL  PERIOD:' 
READ<5»*)TN 

WRIT£<6.*>'  ENTER  PEAK  RESISTANCE*.' 
READ(5.*)RMX 
READ<5» 10) IDUMMY 
10  FORMAT (A2) 

RETURN 

END 


SUBROUTINE  TITLOP< TITLE ) 

IMPLICIT  I NT£GER*2 < I -N » ♦ > 

INTEGERS  T1TLE<72> 

URITE<6»*> '  ENTER  NEW  TITLE:  (72  CHAR.  MAX.)' 
REAB<5»10HTITLE<I)  .1*1.72) 

10  FORMAT (72A1) 

RETURN 

END 
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APPENDIX  B 

EXAMPLE  COMPUTER  SOLUTIONS 
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ft 


«kun  anal/anal  txq 

NORMALIZED  RESPONSE  OF  AN  SOOF  ELASTIC-  PLASTIC  MOOEL 
NOTl  OEFAULT  UNITS  ARE  INCHES.  POUNOS.  AM)  SEC 

unless  otherwise  specified 


ERVISE  SPECIFIED 


available  OPTIONS 

X  -  EXIT  PROCRAM 
A  -  ANALYZE 
0  -  01 SPLAY  OATA 
H  -  HELP 

i  •  initial  pass 

L  -  LOAD INC  FUNCTION 

N  •  numerical  parameters 

0  -  OUTPUT  DEVICE 
P  -  PLOT 

s  •  set  values  for  analysis 
t  -  title 

^SELECT  OPTION.  IX  ,A  J)#,!  A  ,N  A  ,S  ,T  I 

ENTER  NEW  TITLE.  172  CHAR.  MAX.) 

SOLUTION  FOR  B.S  SEC  STRUCTURE  AT  20  KT ,  2S  PSl 
SELECT  OPTION.  IX  ,A  fi  ,H  ,1  A.  M  JO  f  .S  ,T  I 

‘SELECT  SPEIOCR-BROOE,  TRIANGULAR.  OR  EXIT.  IS.T.XI 

ENTER  WEAPONIRT I 
20 

ENTER  PEAR  OVERPRESSURE. 

2S 

ENTER  HEICHT  OF  8URSTIRFEET 1 : 

0 

SOLUTION  FOR  0.S  SEC  STRUCTURE  AT  20  XT  *  25  PS! 

NUMERICAL  PARAMETERS: 

BETA  -  0  10887 

INTEGRATION  INTERVAL  -  0.10000E-01  OF  PERIOD 
INTEGRATION  INTERVAL  -  0  10080E-01  OF  DURATION 
TOLERANCE  -  0. 10000E-02 
NUMBER  OF  PRESSURE  STATIONS  -  2S0 

LOAO  DESCRIPTION. 

WEAPON  IRT !•  20.000 

PEAR  OVERPRESSURE  -  2S.000 

POSITIVE  PHASE  OURATION  -  0.44220 

TIME  OF  ARRIVAL  •  0.S1380 

RANCE  IRFEETI  -  1.7741 

HEIGHT  OF  BURST  IRFEET I  -  0.00000 

PRESSURE-TIME  HISTORY?  CY/Nt 


0H080O 

-3.82760 

13  11820 

2.71724 

0;  10000 

5,45806 

13.14083 

1  70216 

0;- ii  iB0 

*3750045 

13  15424 

0  88924 

0.1  1400 

3  72202 

13.15872 

0.00603 

.0  1 1700 

3.85350 

13  15448 

-085453 

0  12000 

3.08508 

13  M17S 

-1 .69551 

0. ‘12300 

4  11640 

13.12067 

-2  51780 

01 12600 

4.24747 

13  00147 

-3.32178 

0  12000 

4  37820 

13  05432 

-4  10670 

0. 13200 

4.S08S3 

13  00042 

-4  87475 

0  13S00 

4.63837 

12.05601 

-5  62636 

0  13000 

4  76764 

12  80608 

-6  36038 

0  Ml 00 

4.80828 

12.82078 

-7  07901 

0  14400 

5  02421 

12  75547 

-7  78268 

0  H700 

5  15137 

12  67421 

-8  47061 

0.15000 

5.27768 

12  58613 

-0  14413 

0  15300 

5.40307 

12  40130 

-9  80402 

0  15600 

S  52748 

12  30912 

-10  44078 

0  15000 

S.6508S 

12  28246 

-11  08204 

0  16200 

5  77311 

12  16854 

-It  70179 

0  16500 

5.60420 

12  04840 

-12  30886 

0. 16800 

6. 01406 

11  02243 

-12  00323 

0  17100 

6.13263 

11  70048 

-13.4861 t 

0  17400 

6.24085 

11  65277 

-M  05756 

9.17700 

6.36567 

11  50013 

-M  61708 

0  18000 

6.46002 

11.30048 

•IS  16507 

0. 18300 

6  50286 

1 1.20612 

-IS  70458 

0  16000 

6.70412 

11  04044 

-16  23102 

0  18000 

6.61377 

10  88154 

-16  74044 

0  10200 

6.02174 

10.71150 

-17  25730 

0  10500 

7.02708 

10.53644 

-17  7SS2S 

0  10800 

7  13245 

10  3S641 

-18.24373 

0  20100 

7  23500 

10  17161 

-18  72333 

0  20400 

7  33587 

9  08202 

-10.10380 

8.29700 

7.43472 

0  78778 

•19.65533 

0  21000 

7.S3161 

9.58896 

-20.10864 

0  21300 

7  02618 

9.3856S 

-20.SS363 

0  21600 

7.71031 

9  17703 

-20.90013 

0.21000 

7  81003 

8  96588 

-21  41808 

0.22200 

7.80801 

0  74959 

-21  84026 

0  22500 

7.08501 

8.52912 

-22  25343 

0.22800 

8  06018 

8  30455 

-22  65948 

0.23100 

8.15108 

8  07596 

-23.85852 

0  23400 

8.23068 

7  84342 

-23 . 45009 

0  23700 

8.30704 

7.60609 

-23  83485 

0.24000 

8.38281 

7.38676 

-24  21307 

0  24300 

8.4S526 

7  12277 

-24  58446 

0.24600 

8.52525 

6.87510 

-24.94934 

0  24000 

8.S0275 

6  82381 

-25  30811 

0.25200 

8.65772 

6.30807 

-25.66061 

0 . 25500 

8.72012 

6.11003 

-26  00687 

0.25800 

8.77002 

5.84886 

-28.34740 

0.26100 

8.83700 

5.58371 

-26  68219 

0.26400 

8  80158 

5.31524 

-27.01100 

0 . 26700 

8.04338 

5  04352 

-27.33444 

0.27000 

8.00244 

4 . 76858 

-27.65260 

8.27300 

0.03874 

4.40040 

-27.06504 

f. 

t 
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1S0  200  250  300  350  400  450  500 

TIME  (MS) 
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NOTATION 


f  Forcing  function/maxiaua  force 

|  F  Forcing  function 

I'  Fj  Maximum  force 

k  Elastic  function 

|  m  Structural  mass 

P  Peak  overpressure 
s 

r  Structural  resistance 

I  r  Maximum  resistance 

m 

t  Time 

Positive  phase  duration  of  load 
t.  t  Time  to  maximum  response 

ID 

T  Natural  period 
x  Displacement 
B  x  Velocity 

x  Acceleration 
xm  Maximum  displacement 

B  x  Yield  or  peak  elastic  displacement 

p  Newmark  p  numerical  parameter 
At  Time  increment 

-  0  Normalized  displacement 

i"!  Normalized  velocity 

Normalized  acceleration 
*  M  Ductility 

£  Normalized  time  parameter 
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